%将采样集xcell,vcell在谱空间中处理，输出简正频率和模式上的能量
function [energy, omega] = spectrum_analysis(xcell, vcell, tau, en0)
    fs_x = Fourier_Series(xcell, tau);
    fs_v = Fourier_Series(vcell, tau);
    
    sam_len = length(xcell{1});
    ome = (0:sam_len)*(2*pi)/sam_len;
    
    fs_pot = zeros(1,sam_len);
    for j = 1:sam_len
        fs_pot(j) = 0.5 * (ome(j)^2) * fs_x(j);
    end
    fs_kin = 0.5 * fs_v;
    
    [peak_pot, loc_pot] = find(fs_pot, 'minpeakheight', en0*0.01); %寻找势能峰值超过0.01*en0的峰值
    [peak_kin, loc_kin] = find(fs_kin, 'minpeakheight', en0*0.01); %寻找动能峰值超过0.01*en0的峰值
    
    energy = zeros(1,sam_len);
    omega = zeros(1,sam_len);
    
    for k = 1:sam_len
        energy(k) = peak_pot(k) + peak_kin(k);
        omega(k) = 0.5 * (ome(loc_pot(k)) + ome(loc_kin(k)));
    end

    return
end